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Abstract 

Recent studies of the two-dimensional, classical XY magnet in a magnetic field suggest that 
it has three distinct vortex phases: a linearly confined phase, a logarithmically confined phase, 
and a free vortex phase. In this work we study spin-spin correlation functions in this model by 
analytical analysis and numerical simulations to search for signatures of the various phases. In 
all three phases, the order parameter is nonzero and < cos(0(ri)) cos(^(r2)) > remains nonzero 
for r = |ri — oo, indicating the expected long range order. The correlation function for 

transverse fluctuations of the spins, C(r) =< sin(0(ri)) sin(0(r2)) >, falls exponentially in all three 
phases. A renormalization group analysis suggests that the logarithmically confined phase should 
have a spatially anisotropic correlation length. In addition, there is a generic anisotropy in the 
prefactor which is always present. We find that this prefactor anisotropy becomes rather strong in 
the presence of a magnetic field, masking the effects of any anisotropy in the correlation length in 
the simulations. 

PACS numbers: 64.60.-1, 64.60.Cn, 75.10.Hk 
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I. INTRODUCTION 



The two-dimensional classical XY model is important because it serves as a paradigm for 
many systems, including easy-plane magnets, two dimensional solids, thin film superconduc- 
tors and superfiuids, magnetic bubble arrays, and certain one-dimensional quantum systems. 
Of particular interest in these systems is the behavior of their topological excitations, which 
in the XY model are vortices. In spite of the Mermin-Wagner-Hohenberg theorem P], which 
forbids long-range order (i.e., a non- vanishing order parameter) at any finite temperature 
in two dimensions, a (Kosterlitz-Thouless, or KT) phase transition occurs in this system, in 
which the vortices go from a bound pair phase at low temperature to an unbound vortex 
phase at high temperature j^. While the order parameter vanishes in both phases, corre- 
lation functions show different behaviors in the different phases. In the high temperature 
3hase, the correlations decrease with increasing distance according to an exponential law 
^1 El 'while in the low temperature phase, the correlations decrease according to a power 




A non- vanishing order parameter can trivially be restored by the application of a magnetic 
field which tends to order the XY spins. The behavior of the vortices in this situation 
has recently been revisited via renormalization group (RG) and simulation studies, 

leading to the conclusion that vortices can still unbind due to thermal fluctuations in this 
system, in a two-step process. At the lowest temperatures, vortex-antivortex pairs are 
linearly confined by a string of overturned spins. As temperature is increased, these strings 
first undergo a proliferation transition, but the vortices remain confined due to a residual 
logarithmic attraction. With increasing temperature this attraction is overcome and the 
vortices deconfine. Importantly, unlike the KT transition, these are not phase transitions 
in the thermodynamic sense, as they do not lead to singularities in the free energy 0|. 
Nevertheless, the transitions can int^duca .„gulant,e= .n co.elat.on functions i, lead 
to dramatically different behaviors for the physical systems in which the XY model with a 
magnetic field is realized. This includes the behavior of bilayer thin film superconductors , 
the bilayer quantum Hall system ^ , and bosons in a linear optical lattice, tunnel-coupled 
to a bulk superfiuid reservoir 

Since a standard way to characterize phases and phase transitions is via the behavior of 
correlation functions, it is natural to search them for signatures of the deconfinement transi- 
tions. In this paper, we study spin-spin correlation functions of the system in all three phases 
through both analytical analysis and numerical simulations. Due to the external field, all 
the three phases will have a non- zero order parameter, and for a symmetry-breaking term 
(i.e., magnetic field coupling) of the form — /icos^(r), with 6{r) the angular variable of 
the XY spin located at lattice site r, one easily confirms that the correlation function 
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< cos(6'(ri)) cos(^(r2)) const., for r = |ri — oo, indicating long range order. A 

more interesting correlation function C(r) =< sin(6'(ri)) sin(6'(r2)) > measures the fluctu- 
ations of the spins transverse to the direction of the symmetry-breaking field. We find in 
all three phases C(r) ~ e~^^'^, but that in the log-confined (Log) phase, A = X{a), with 
a the angle between r and the symmetry axis of the underlying lattice. (Throughout this 
study, we will assume the XY spins reside on a square lattice.) This anisotropy in principle 
is a signature that distinquishes the Log phase from the linearly confined and deconfined 
phases. We will see, however, that there is in addition an anisotropy that exists in all three 
phases (and in the deconfined phase of the XY model without a symmetry-breaking field) 
which renders the identification of the anisotropy due to the vortex phase very challenging 
in numerical simulation. 

The organization of this article is as follows. In the next section, we study correlation 
functions in the three phases via some simple analytical models. Section III reports the 
results of our numerical simulations. Section IV discusses the anisotropy found in our 
simulations for all the phases, which we show is a ubiquitous effect of the underlying lattice. 
In Section V we report on attempts to distinquish the lattice anisotropy from that associated 
with the vortex phase in the numerical data. Conclusions are presented in Section VI. 



II. CORRELATION FUNCTIONS: ANALYTICAL MODELS 

The Hamiltonian for the XY model with a magnetic field is 

H = -i^ E<r,r'>cos[^^(r) - ^(r')] - /iE.cos[^(r)], (1) 

where < r, r' > refers to nearest neighbor sites on a two dimensional square lattice. As 
discussed in the Introduction, we focus on the order parameter and transverse fluctuation 
correlation functions < cos[6'(ri)] cos[^(r2)] > and < sin[6'(ri)] sin[0(r2)] >. Both these can 
be computed if we can calculate < >. A generating functional that allows us to 

do this is 

2 

= J DOe-^+'^r J{r)0{r) I J Dee-"", 

where J(r) = 5{v) gives the order parameter average, and J(r) = 5{y — vi) ± 5{y — V2) 
generates < e*^('"i)=^*^(''2) >. 

We may study our problem in the Villain model This essentially involves 

replacing e'^^°^^ with X]m=-oo ^"'^^^"^'^'"^^''^ wherever it appears in the partition function. 
The main idea is that the long distance physics depends only on the symmetry of the 
Hamiltonian, not the details of the interactions. 
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Following this prescription we may write 

= Es , T.Tr /[L'^]e"27r^^-'"5k^'''e^E(v-s+T+j)e{r)_ 
Integrating out 9, we have 

Z[J]=Es^^,Ene~^^^'^''^^^'s{^-S + T + J). (4) 

Note that the (5-function appearing here is actually a Kronicker delta with integer arguments. 
The delta function allows a further simplification if we write 

S^, = e^^d^n + + rj^, 

V ■ A= T (5) 
V-r] = J, 

where n, A, rj are integer fields. (Care must be taken in choosing an allowed form of A such 
that all configurations of the integer field S are correctly produced [^.) From Eqs. (4) and 
(5), we arrive at 

1 Y^l. a„iA i„|2 1 



Z\J] = y . e"27r^'^-"''"+^'^+'"'l'-27^^l^-^^ (6) 

L J ^ — '71, All, 



A. Linearly Confined Pliase 

In the linearly confined phase, the integerness of n and A is unimportant [6|, and we may 
treat them as continuous fields. Integrating them out, we find the generating functional 
takes the form 

where ^ = \jK/h. i has the interpretation of the width of a string connecting a vortex- 
anti-vortex pair 0,0. 

• Order parameter 

For the order parameter, J(r) = 5(r), \J{q) \ = 1, and we get 

• Correlation function 
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(10) 



For < e*^*^""!^ iS{r2) y, take J(r) = S{r — ri) — S{r — r2), and using the fact that 
< sin(6'(ri)) cos(^(r2) >= by symmetry, we find 

^ 1 - cos(gr) 

< cos(6'(ri)) cos(6'(r2)) - sin(6'(ri)) sin(0(r2)) >= e 

where r = |ri — r2|. 

Similarly, for < e«^(ri)+««'(r2) ^ ^g^j^g _ ^(^j. _ j.^^ _|_ ^(^j. _ j.^-)^ g^^^ 

^ 1 + cos(gr) 

< cos(6'(ri)) cos(6'(r2)) + sin(6'(ri)) sin(6'(r2)) >= e 
Thus we have 

1 . ,2 1 - cos(gr) 1 . ,2 1 + cos(gr) 

<cos(^(ri))cos(^(r2))>=^[e + 9 + e + 9 ] ^^^^ 

and 

1 2 l-cos(gr) 1 2 l + cos(gr) 

C(r) =< sin(^(ri))sin(^(r2)) >= ^[e 1/4 _e ]. 

(12) 

For large distance r, we have approximately 

cos(gr) 



< cos(e(rO) cos(^(r2)) >^ I - ^ J d^^pf^ 
C{r) ^< sm{e{n))sm{e{r,)) >^ ^ I d'q^^^fl, 



(13) 



C{r) ^ Ko{r/0 ~ ^p^'/?, (14) 

where Kq is the 0-th modified Bessel function, and the last form of C{r) is valid at large 
distances. 



B. Logarithmically Confined Phase 



In this phase, a suitable Hamiltonian for the long distance physics which preserves the 
underlying lattice symmetry ^ is 

H=^I dMbl + Kpql]\<P,\' + [ql + Kpql]\ct>,\^ + e^|g.g,(0i + 02)^}, (15) 

where the ii"pg^|0ip, -ft'pg^|02p terms are generated by the renormalization group fiows of 
operators that model the integerness of the fields from which the continuous fields 0i,02 
were derived. The generating function < e^^ ■^^ > can be calculated by the path integral 

Z[J] = f D<P,D<P2 exp{-27^ / d\{[ql + Kpq'y]\<P,\' + [g^ + Kpql]]^^^]' + eM<Pi + + J?}}- 

(16) 
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Integrating out the fields 0i, 02, we have 



Z[J] = exp\-4rrId^q ^J^^k — 4 — T^\, , , 

?y 4\, 22 



SO that for large r, 



^/ X r ,2 cos(gr) 

2A' J " Kp(si + 4-i,i^ - (18) 

We see that in the Log phase, the mass terms oc Kp generated by the renormalization group 
flows introduce a characteristic anisotropy in the correlation function C(r). If evaluated 
at large distance, C(r) depends exponentially with distance, but the correlation length is 
angle-dependent due to the Kp term. 



C. Free vortex phase 

The free vortex phase occurs at high temperature, where the coefficients K and h are 
small. The simplest way to understand the form of the correlation functions is to return to 
the orginal XY model and compute the correlations functions perturbatively. The transverse 
spin fluctuation correlation function is defined as 

C(ri -rs) =< sin(^(ri))sin(^(r2)) >= J[D9] sin(^i) sin(^2)e^^<'-.'-'> . 

(19) 

To second order in K and h we may expand this as 
< sin(^(ri)) sin(e(r2)) >= J[De] sin(ei) sin(02) 

X n<r,r'>{l + Kcos{e{r) - ^(r')) + ^[cos{e{r) - e{r'))]^} JJA^ + /icos^(r) + ^ cos^ ^^(r)} 

(20) 

where Pj refers to self- avoiding paths connecting ri and r2, and \Pj \ is the length of the path. 
The shortest path gives the dominate contribution, so that C{r) ~ e~^^^, with correlation 
length A = ^ {2/K) — h'^/8 " '^^^^ same result can be obtained by approximating the sum 
over self-avoiding paths as a sum over all random walks, — > Yl'p ■ The summation over 
random walks can be calculated by a generating function method [l^. We first note that 
C(r) can be written as 

^ ^ (21) 

= E^^iv(r)x^ 
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where x = K/2 + Kh"^ /16, r = r2 — ri = {m,n), {m,n) are the coordinates of r and 
P/v is the number of path of N steps. P/v may be conveniently written in terms of the 
function f{zi, Z2) = zi + z^^ + Z2 + Z2^ , which obeys the identity (/)^ = Ylik i ^k,i^i^2i where 
PAr(r) = PN{m,n) = c^^. Thus we have 

PN{r) = <„ = -^Jdz,dz2zT-'zr\f{z)r. (22) 
Using eq.(21) and eq. (22) and setting zi = e*^^, Z2 = e^'^y , we obtain 

C(r) ^ Pnx"" = J yt^i o I T'L r^- (23) 

V / J (^2ny 1 - 2x[cos(g^) + cos(gy)J ^ 

For small x and large r, the integral may be evaluated leading again to the result C{r) ~ 

Using similar methods, the order parameter and order parameter correlation function 
may be shown to behave as 

< cos(6') >~ ^, 

< cos6'(ri) cos6'(r2) >~ ^ + e"''/^, (24) 

A = 1 

ln{2/K) -3h''/8' 

Finally, we note that similar results may be obtained within the Villain model, although 
the derivation is considerably more involved than this. 



III. NUMERICAL SIMULATION 

In the previous section, we found a characteristic anisotropy in the Log phase that is not 
present in either the linearly confined or unbound vortex phases by computing correlation 
functions near fixed points of these phases 0, 0]. We now search for this behavior in 
simulations. We shall see that anisotropy in the correlation functions is quite generic in the 
XY model, presumably due to irrelevant operators that were not included in our effective 
Hamiltonians. Our results suggest that very high precision simulations are needed to detect 
the effects of the vortex phase in the correlation functions. 

We employ a standard Langevin dynamics approach for our simulations. Angular vari- 
ables < 9 < 2n are represented on a square lattice with periodic boundary conditions. 
The model is equivalent to that of an array of superconducting grains, treated within the 
resistively shunted Josephson junction model 0|. The equations of motion for the angular 
variables are 

rf^g(r) _ SHxY I de{r) . 
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The spin stiffness was fallen as = 1, setting our unit of energy. The effective moment of 
inertial for the spins F were also taken to be 1, setting the unit of time in our simulation. 
The viscosity t] was taken to be 0.1, which we found yields good equilibration times and 
statistics. C(r, t) is a random torque satisfying < C(r, t)C(r', t') >= 2riT6j.y6(t — t') with T 
being the temperature of the system (chosen to be 1.2 for the results discussed below.) The 
system size is = Ny = 61. Our Hamiltonian Hxy takes the form 

HxY = -KE<r,r'>cos[e{r) - e{r')] -hj:^cos[eir)], (26) 

with K = 1. A typical run consists of 5 x 10"^ time steps. Each time step is 0.1 (in units 
of ^yV/K). We also eliminated the initial 2 x 10^ steps in these runs for equilibration. We 
repeated runs for each set of parameters with several different seeds, allowing us to estimate 
the statistical error. 

Our simulation results are summarized in Fig. 1-4. The curve with solid dots (squares) 
is for the transverse spin fluctuation correlation function C(r) in the x (y) direction. Stars 
show correlation functions along a direction with angle ir/A with respect to the x axis. Fig. 
1 shows results for h=0, while Figs. 2, 3 and 4 are respectively for h = 0.417, 0.25 and 0.15. 
These values were chosen to place the system in the linearly confined. Log, and unbound 
vortex phase, respectively, as found in an earlier study in which the vortex phase was probed 
by measuring fluctuations in the vortex dipole moment of the system 3|- 

In all the simulations, we see clearly that C(r) ~ e^^'^'^,r = |ri — r2|; i.e., the correlation 
function falls exponentially at large distances. This agrees with our expectations from 
Section II. The correlation functions in the x and y directions are the same within the 
numerical accuracy of the simulations, as expected from the lattice symmetry. There is 
also a small difference between the amplitudes of the correlation functions along the bond 
direction and along the diagonal. However, the form of this anisotropy is not what we expect 
from Section II, which suggested anisotropy in the correlation length in the Log phase. This 
should lead to different slopes for C(r) in the log- log plots along the bonds and the diagonals 
at large distances. Although there is some evidence for this, it occurs only at very large 
distances, where our statistics are relatively poor and an accurate measurement was not 
possible. 

That most of the anisotropy appears as an offset in the correlation functions in the log-log 
plots suggests that our results are dominated by the prefactor of the exponential function 
for C(r). We explain the origin of this in the next section, and argue that it is a generic 
behavior that is independent of the vortex phase of the system. 
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IV. LATTICE EFFECTS 



To take the lattice effects into account, the integration in eq. (7) 
should be modified to 

f (no) 

^l/e + A-2cos{q,)-2cos{qy) ^ ^ 

where our unit of length is the bond length. Following the methods of Section II, the 
resulting transverse spin fluctuation correlation function becomes 

^('■' ~ ^ '''' v^rraffc:^- (29) 



For large distances, C(r) may be calculated using the stationary phase approximation 
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17l |. with the result 



expl-q^TTi:, - qyTUyi 



Here r = (m^,, m^), and g*, q* are determined by 

1/^2 _ 4 - 2 cosh(g*) - 2 cosh(g;) 
mj^sinh(g*) = ma;Sinh(g*), 

and R is given by 



(31) 



^ _ y sinh (g*) + sinh (g^)) ^^2) 

I sinh^(g;) cosh(g*) + sinh^(g*) cosh(g;)| ' 



For large ^, the solution to Eq. EH is 

qx = I cos a 
qy = ^ sin a, 

where tana = rriy/mx- In this case, the correlation length is A = ^ and 



(33) 



^ l + (a/0'sin2(a)cos2(a)' ^^^^ 

Thus we see that the lattice structure induces anisotropy in the prefactor of the correlation 
function. The ratio of the correlation function along two different directions, with angles 
«!, a2 with respect to the x axis, has the form 
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TABLE I: Fitting parameter I 



h 


A 




Kp 


a 


h 


0.0 


0.5165 


4.6173 


-0.1940 


-0.8616 


0.3821 


0.25 


0.08154 


1.7834 


-0.03805 


-0.1673 


0.3765 


0.417 


0.05268 


1.4154 


-0.007307 


-0.1095 


0.3060 



for a/^ << 1. This prefactor anisotropy persists even at very long distance, and shows up 
as an offset in a plot of C(r) along different directions when displayed in a log-log plot. This 
is exactly what we see in our numerical simulation, and it is apparent that for the distances 
where we can make reliable measurements, it is the dominent effect. We note that, although 
this result was presented for the linearly confined phase, one may show that in all the phases 
discussed in Section II, a lattice regularization of the momentum integrals leads to similar 
effects. 



V. DISCUSSION 

As discussed above, there is anisotropy in the correlation functions for all three phases 
due to lattice regularization. One can search for extra anisotropic behavior in the Log phase, 
as in the section II. To separate the two kinds of anisotropics, we fit our numerical results 
by the formula 

C{r)^Ajd^q . , , 

,,.22 Kp^i^y , 4 , , 4 4 36 

q " 

We choose the general form ^— ^ (compatible with lattice symmetry) to mimic the 

anisotropic behavior apart that should not arise due to the lattice effects but is expected 
in the Log phase. The term aq^ is a higher order isotropic term and the term bq^qy is 
the anisotropic term due to lattice effects, expected from expanding the cosines in a lattice 
regularization, and which for large enough distances should be the leading order corrections 
due to this regularization. We fit our numerical curves for the correlation function within a 
suitable distance interval (generally from 3 to 6 correlation lengths.) If the distance is too 
small, short distance correlations become important and presumably need further terms to 
be correctly modeled. On the other hand, when r is too large, C(r) becomes very small, 
and numerical errors spoil the fit. The resulting fitting parameters are shown in Table I. 

We see that Kp is small but nonzero for all three phases, and moreover is largest when 
h — when we expect no such term. We believe this is due to the fact that one needs to 
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TABLE II: Fitting parameter II 



h 


A 




a 


h 


0.0 


0.5163 


4.6207 


-1.1499 


3.0605 


0.25 


0.08097 


1.7884 


-0.1800 


0.5610 


0.417 


0.05256 


1.4167 


-0.1110 


0.3311 



measure the correlation function at much larger distances than we can achieve numerically 
to see its true asymptotic behavior. This is consistent with earlier simulations 0] which 
required detection of vortex-antivortex pairs at effective separations of order the system 
size, typically ~5 times larger than the maximum r used in our fits. To support this point 
of view, we fitted our numerical results by Eq. EHl but set Kp = by hand; i.e., 

C{r) = AJ rf^g— (37) 
1 + ^ q +aq + bq^q^ ^ ' 

The obtained fitting parameters are listed in Table II. 

To further compare the fits with one another and the numerical data, we use Eqs. EHl 
and|nZlto calculate model correlation functions over long distances with the fit parameters. 
The difference between these turns out to be extremely small even at large distances. (For 
h = 0.0, it is 0.81%; for h=0.25, it is 4.9%; for h=0A17, it is 4.0%.) The absolute difference 
is 3.3 X 10"^ for h =0.0; 2.8 x 10"^ for h=0.25; and 4.9 x lO'^ for /i=0.417. Most importantly, 
in all cases these differences are of the same order or smaller than our numerical error. To 
demonstrate this, in Fig. 5 we compare the correlation functions C(r) (both along a bond 
direction and along the diagonal) for h=0.25, obtained by numerical simulation, from Eq. 
EHl and from Eq. EZI The curves from the three methods are nearly indistinquishable. In 
the inset, we show the correlation function for an arbitrary direction, tan(a) = 1/2, where 
a is the angle with respect to the x axis. Again, we can hardly see any difference among 
the curves. Clearly, much higher precision data than we were able to generate is required to 
distinquish the effects of lattice anisotropy from the effects of logarithmically bound vortices. 



VI. CONCLUSION 



We have discussed the effect of the various vortex phases of the XY model in a magnetic 
field on spin correlation functions, focusing on transverse spin fluctuations. Calculations 
based on the fixed point Hamiltonians suggest that the Log phase should have a unique 
signature in an anisotropic correlation length, but there is anisotropy common to all the 
phases due the lattice realization of the model, that shows up as irrelevant operators in the 
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renormalization group. Distinquishing these lattice effects from those of the vortex phase 
turns out to be quite difficult for system sizes and run times that are readily attainable 
in our Langevin dynamics simulations. While we believe the anisotropy is present, direct 
confirmation will require higher accuracies and larger systems than were achieved in our 
own study. 
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FIG. 1: Correlation functions vs. distance in a log-linear plot for ^=0, system size=61. Averaged 
over 9 seeds. Curves with solid dots and squares are for correlation functions in the x and y 
directions, and cannot be distinquised in the figure, indicating the errors are smaller than the 
symbol sizes. Curve with stars is for correlation function in the direction with angle 7r/4 with 
respec to x axis. 
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FIG. 2: Correlation functions vs. distance in a log-linear plot. h=0A17 (linearly confined phase), 
system size=61. Averaged over 7 seeds. Curves with solid dots and squares are for correlation 
functions in the x and y directions, and cannot be distinquised in the figure, indicating the errors 
are smaller than the symbol sizes. Curve with stars is for correlation function in the direction with 
angle 7r/4 with respec to x axis. 
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FIG. 3: Same as Fig. 2, except h=0.25 (Log confined phase), and data is averaged over 5 seeds. 
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FIG. 4: Same as Fig. 2, except ^=0.15 (free vortex phase), and data is averaged over 2 seeds. 
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FIG. 5: Correlation functions vs. distance in a log-linear plot, h =0.25 (Log confined phase), 
system size=61. The two linesin the main panel (for x direction and the direction with angle 7r/4 
with respect to x axis are actually two pairs of three overlapping curves: the results obtained from 
numerical simulation, fitting formula Eq. 1311 and fitting formula Eq. 1371 demonstrating that our 
data do not distinquish the expected anisotropy in the Log phase from lattice anisotropy effects. In 
the insert are correlation functions in the direction with angle tana = 1/2. Two curves obtained 
from Eqs. OH and ISZI are again indistinquishable. 
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